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Aims. We present simulations of the perfomances of the future SPHERE IFS instrument designed for imaging extrasolar planets in 
the near infrared (Y, J, and H bands). 

Methods. We used the IDL package code for adaptive optics simulation (CAOS) to prepare a series of input point spread functions 
(PSF). These feed an IDL tool (CSP) that we designed to simulate the datacube resulting from the SPHERE IFS. We performed 
simulations under different conditions to evaluate the contrast that IFS will be able to reach and to verify the impact of physical 
propagation within the limits of the near field of the aperture approximation (i.e. Fresnel propagation). We then performed a series of 
simulations containing planet images to test the capability of our instrument to correctly classify the found objects. To this purpose 
we developed a separated IDL tool. 

Results. We found that using the SPHERE IFS instrument and appropriate analysis techniques, such as multiple spectral differential 
imaging (MDI), spectral deconvolution (SD), and angular differential imaging (ADI), we should be able to image companion objects 
down to a luminosity contrast of ~ 1CT 7 with respect to the central star in favorable cases. Spectral deconvolution resulted in the most 
effective method for reducing the speckle noise. We were then able to find most of the simulated planets (more than 90% with the 
Y-J-mode and more than the 95% with the Y-H-mode) for contrasts down to 3 x 10~ 7 and separations between 0.3 and 1.0 arcsec. 
The spectral classification is accurate but seems to be more precise for late T-type spectra than for earlier spectral types. A possible 
degeneracy between early L-type companion objects and field objects (flat spectra) is highlighted. The spectral classification seems 
to work better using the Y-H-mode than with the Y-J-mode. 

Key words. Instrumentatiomspectrographs - Methods: data analysis - Techniques: imaging spectroscopy - Stars: planetary systems 



1. Introduction 

A large number of extrasolar planets have been discovered in 
the last fifteen years through indirect methods such as radial 
t-H ■ velocities and transits . Although in the past few years some 
,— 1 , objects with planetary mass have been imaged around stel - 
es! ■ lar and sub stellar objects like HR 8799 dMarois et al. I [2008). 
• i-h ; Fomalhaut dKalas et al. Il2008h. 2 M 1207 dChauvin et aLT 2009). 
. and p Pictoris (iLagrange et al. 1 1201 Oh . imaging of extrasolar 
H ' planets is still very challenging because of the high planet vs 
. . .' star luminosity contrast (10~ 6 for young giant planets and down 
to 10~ 8 -10~ 10 for old giant and rocky planets) and the small sep- 
aration with respect to the central star (few tenths of arcsec for a 
planet at ~10 AU at some tens of pc). 

The next generation of instruments aimed at imaging extrasolar 
planets will exploit extreme adaptive optics (XAO) systems to 
correct aberrations up to a high order, providing a high Strehl 
ratio (SR) and high-efficiency coronagraphs to attenuate the on- 
axis PSF and reduce its diffraction pattern. The combination of 
these two devices should be able to reduce the stellar background 
down to a value of around 10~ 5 at separations of a few tenths 
of arcsec. The residual background will be given mainly by the 
speckle noise generated by the atmosphere and the telescope 
pupil-phase distortion. To further improve the contrast achiev- 



able with these instruments, it will be mandatory to apply differ- 
ential ima ging techniques, s uch as angular differential imaging 
(ADI) dMarois et al. 1 12006). sim ultaneous sp ectral differential 
imaging (S-SDI) (see e.g. | Marois et al. 1 20051) . and spectral de- 
convolution (SD) (see lThatte et al. 1120071) . ~ 
In the next years in particular, three instruments will be able 
to exploit these techniques to image extrasolar planets. These 
are the G emini Planet Imager ( GPI) at the Gemini South 
Telescope (iMacintosh et al. II2006I). SPHER E at the ESO Very 
Large Telescope (VLT) dBeuzit et al. 1 12006). and Project 1640, 
which is already working at the 5 m Palomar telescope (see 
ICreppet al. Il2010h . 

In particular, SPHERE will include three scientific channels: (i) 
a differential imager and dual band polarimeter called IRDIS 
that will operate in the n ear infrared between the band Y and 
K s dDohlen et al. 1 120081) : (ii) a polarimeter called ZIMPOL 
that will perform differential imaging exploiting the polar- 
ized light reflected from the plan etary atmosphere in the vi- 
sual band ( Thalm ann et al. 1 12008): (iii) an integral field spec- 
trograph (IFS) that will supply simultaneous images at differ- 
ent wa velengths in the ne ar infrared between the Y and the H 
bands dClaudi et al. II2008I) . 

Integral field spectrographs also have the potential of providing 
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the spectra of the detected faint companions at close separation, 
thus allowing much better characterization. IFSs similar to the 
one designed for SPHERE are also present in GPI and in Project 
1640, and are foreseen for future pla net imagers like EP ICS de- 
signed to work for the future E-ELT (Kaspe r at al. 1120 1 Ob - In this 
paper we present the results obtained from the simulations we 
developed to evaluate and to optimize the performances of the 
SPHERE IFS. In Section[2]we give a very short summary of the 
SPHERE IFS instrument, in Section [3] we describe the simula- 
tion tools that we used in our work, in Section|4]we describe the 
methods used for the data analysis of the output of our simula- 
tions, in Section [5] we present the results of the various simula- 
tion runs, in Section [6] we describe the software that we wrote 
for the IFS data analysis and the results obtained testing it on 
the output of our simulations, while in Section Q we report our 
conclusions. 

2. SPHERE IFS description 

The SPHERE IFS is designed to work in two different wave- 
length ranges: (i) 0.95-1.35 //m (Y-J-mode) with a resolution of 
R=50 and (ii) 0.95-1.65 yum (Y-H-mode) with a resolution of 
R=30. These two ranges and resolutions are ach ieved through 
two different dispersers (two Amici prisms - see IOlival l2000). 
The IFS is composed of several subsystems: 

- the integral field unit (IFU) 

- the collimator optics system 

- a filter wheel 

- the disperser optics system 

- a camera optics system that can be moved to focus spectra on 
the detector or to produce dithering to reduce noise related 
to the flat fielding 

- a 2048x2048 Hawaii II detector with pixel of 18 fim housed 
in a cryostat 

The novel lenslet IFU concept upon which this spectrograph 
is based (BIGRE, lAntichi et al. 1 12009) allows the entrance slits 
plane to be made of images of the telescope focal plane and 
not of images of the telescope pupil, a s ordered in the clas- 
sical TIGER design dBacon et al. 1995). In this design, each 
lenslet is an afocal system with two powered surfaces. The thick- 
ness of the array is then given by the sum of the focal lengths 
of the lenslets of the two arrays. The main advantage of the 
BIGRE configuration over the TIGER one is that it allows a 
strong reduction of the cross-t a lk bet ween adjacent lenslets as 
demostrated by lAntichi et al . (2009). The microlens array is 
composed of 145 x 145 hexagonal lenslets with a pitch of 161.5 
ym (corresponding to ~0.012 arcsec). Each lenslet is masked 
with a circular aperture with a factor of 0.98 to avoid straylight. 
The full field of view (FOV) of the instrument is a square with 
a side of 1.77 arcsec. The total length of the whole instrument 
from the first surface of the IFU to the detector plane is 1061.89 
mm. A more d etailed description o f the whole instrument can be 
found, e.g., in lClaudi et all (1201 Ol) . 

3. Simulation description 

We exploited two software tools for our simulations: 

- the SPHERE package of the CAOS software 

- the CSP coddj 

1 CSP stands for CHEOPS Simulation Program. The code was orig- 
inally developed in the context of the first feasi bility study for a direct 
imaging planet finder for VLT called CHEOPS dFeldt et al. Il2003h . 



The CAOS system (ICarbillet et al. 1120041) is an IDL based soft- 
ware that aims to simulate the behavior of a generic adaptive op- 
tic (AO) system from the atmospheric propagation of light to the 
sensing of the wavefront aberrations and the corrections through 
a deformable mirror. This is done with a Fraunhofer approach, so 
that it cannot be used to properly evaluate the impact of Fresnel 
propagation (see Section lBTTT i. An end-to-end numerical tool has 
been developed for the simulation of the whole SPHERE in- 
strument within the CAOS environment. It contains detailed in- 
strumental modeling of the Extreme adaptiv e optics systems, of 
IRDIS and ZIMPOL (Carbille t1st1iil2008l) . A module simulat- 
ing the SPHERE IFS has been also developed to properly take 
both the real and the imaginary parts of the image forming on 
the lenslet plane into account. In principle, this could allow a 
complete treatment of the cross-talk among the lenslets when 
studying the impact of light propagation through the BIGRE. 
However, the execution of this module turned out to be very time 
consuming so that it was not possible to use it for a large number 
of detailed simulations. To overcome this difficulty, we used a 
shorter code that calculates the impact of the cross-talk between 
adjacent lenslets (coherent) and adjacent spectra (incoherent) by 
providing the beam propagation over a sub-samp le of 7 hexag- 
onal le nslets. This code is described in detail in lAntichi et al. I 
(120091) . After running this code we concluded that a value of the 
cross-talk equal to or less than 10 -2 was completely adequate for 
meeting the objectives of our instrument. 

We then decided to use our (IDL oriented) code called CSP to 
perform all the simulations of light propagation within the IFS, 
while we decided to use the SPHERE CAOS package to pro- 
vide real intensities over the IFU entrance focal plane as input 
for CSP. For this, we performed simulations using the CAOS 
IRDIS module with 100 atmospheric phase screens at 64 differ- 
ent wavelengths ranging between 0.95 and 1.35 fim (or between 
0.95 and 1.65 ymi in the Y-H-mode case). There are enough at- 
mospheric screens is large enough to ensure that static speckles 
dominate noise, as expected in real cases, and to ensure that the 
PSF has an overall shape representing a realistic stellar halo. 
In Figure Q] we display a monochromatic PSF obtained from 
the CAOS simulations. Even if the SPHERE package of the 
CAOS system allows simulating different types of coronagraphs, 
we p referred to use only a 4-quadrant one for our simula- 
tions dBoccaletti et al. Il2008h . In this way, however, our results 
are still representative as we were not interested in investigating 
the performances of all the SPHERE instrument coronagraphs. 
The choice of making a simulation with only 64 wavelengths 
was determined because for more wavelengths, the program sat- 
urates our computer memory. However, CSP requires 269 PSFs 
at different wavelengths as input and to obtain them, we per- 
formed interpolations starting from the ones resulting from the 
CAOS simulatiotfl 

An ea rly version of the CSP code, described in iBerton et al. I 
(2006), has been deeply modified to take variations in the in- 
strument optical design into account. CSP only considers the 
real part of the image on the lenslet plane and then propagates it 
through the IFS spectrograph using a Fraunhofer approach, but 
it can include a treatment of the cross-talk through a parametric 
approach. The code can be divided into different parts: 



2 Note, however, that 64 wavelengths are enough to properly sample 
both the expected sp ectral extension of s peckles over the whole FOV 
(hypers ampling - see Antichi et al. 2009) and each pixel with approxi- 
mately two sampling per pixel. This should guarantee that interpolation 
errors are under control. 



Mesa et al.: IFS simulations 



3 




Fig. 1. Monochromatic PSF resulting from the CAOS simula- 
tion. The bright corona corresponds to the outer working an- 
gle of the XAO system. Its radius is roughly 0.5 arcsec at the 
working wavelength. This image has been obtained using a 4- 
quadrant coronagraph. The cross structure centered on the center 
of the image is the signature of this type of coronagraph. 



- The image formation part simulates the propagation of the 
light through the instrument and its main goal is to produce 
a final image with all the spectra. For each spectral step, the 
exact number of photons passing through every microlens is 
calculated, as well as the correct position projected on the 
detector of the center of each microimage. The intermediate 
images generated by this process are then convolved with a 
microlens PSF prepared in advance. All the monochromatic 
images are properly shifted to account for the spectral dis- 
persion due to the Amici prisms and are then summed up to 
create the spectra. Finally the code adds noises to the image: 
Poisson noise and all the detector noises. An example of the 
output of this part of the code is given in Figure [2] 

- The calibration part performs the same procedure as is de- 
scribed at the previous point using a monochromatic uniform 
illumination of the IFU (and not a PSF) as input. This part 
simulates the wavelength calibration lamps and is performed 
at three different wavelengths. The code, then, reads from 
the images resulting from these procedures the spectra (hav- 
ing a template that indicates the position of the spectra at 
the minimum wavelength). Every spectrum is fitted with a 
Gaussian curve (using the IDL routine GAUSSFIT) and the 
code finds the center of the Gaussian and its error. Finally, 
through an appropriate interpolation, the code calculates the 
shift for each lenslets, using the positions of the previuosly 
calculated three centers and the theoretical position of the 
center (well known because we know the wavelength of the 
calibration lamp). At the end, these results are saved in a 
wavelength map file where at every pixel of the image is as- 
sociated a well-defined wavelength. 

- The last part of the simulation procedure is to create of the 
datacube with the monochromatic images that will be the 
instrument final output. To this aim we derive a rectangu- 
lar grid from the original hexagonal pattern of the IFU. This 
is done by creating a square grid of points at the same dis- 




Fig. 2. A small portion of spectra resulting from one of the CSP 
simulations. 



tance from each other. For every wavelength, the flux value 
associated to every point of the grid is calculated by con- 
sidering the three nearest points at the given wavelength as 
calculated in the calibration step of the procedure and saved 
in the wavelength map. The calculation is made with a mean 
of the fluxes of these three points weighted according to the 
distance from the grid point considered. 

Where not specified, the simulations were performed assum- 
ing a GO spectral type central star with an absolute magnitude of 
J=3.75 and at a distance of 10 pc from the Sun. A total exposure 
time of 1 hour was generally simulated even if, in some cases, 
longer exposure times were simulated. For these simulations we 
assumed a readout noise of 10 e— , a dark current of 0.1 e— , and 
a flat field error of 10 (hereinafter detector noise). 

4. Data analysis methods 

High-performance coronagraphs within an extreme AO system 
like the one adopted in SPHERE, which gets its sampling fre- 
quency equal to 20 cycles/pupil, allows imaging of compan- 
ion objects down to a contrast of 10 within the whole FOV 
of its IFS (2.5 arcsec diagonal), except for separations smaller 
than ~0.1 arcsec from the central star. However, to fulfill the 
goal of the SPHERE instrument to image giant planets around 
young nearby star, contrasts of ~ 10~ 6 - 10~ 7 are requested. 
To this aim, speckle noise has to be reduced by a further factor 
of 10 - 100. This is done by applying some differential imag- 
ing analysis techniques to the final datacube extracted from the 
IFS data, su ch as the simultane ous spectral differential imag- 
ing (S-SDI) (Maroi s et al. 1120001) and the spectral deconvolution 
(SD) dThatte et al. 1120071) . A natural evolution of the S-SDI tech- 
nique, using an IFS, has been defined during our work and is 
described in detail in Section 14.11 Another possible technique, 
normally applied with ot her analysis techniq ues, is angular dif- 
ferential imaging (ADI) dMarois et al. 1120061) . 
In this section we briefly present the algorithms we used to im- 
plement these methods. 
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4.1. Multiple differential imaging 

As previously said, this techique is an extension to more spectral 
channels of the previous S-SDI techniques, such as the single 
differential imaging for two channels and the do uble differential 
imaging for three channels dMarois et al. Il2000l) . The final result 
of the CSP simulation code consists of a datacube composed of 
33 (for the Y-J-mode) or 38 (for the Y-H-mode) monochromatic 
images. On these images we apply the following steps: 

1. The images are divided into two groups: planetary images 
(monochromatic images at wavelengths where the planet sig- 
nal is potentially present) and reference images (monochro- 
matic images at wavelengths where the planet signal is very 
weak or absent) according to giant-planet atmosphere mod- 
els. 

2. We then distinguish two different cases: 

- Single differences: 

- A reference image is assigned to each planetary im- 
age. 

- For each pair, the reference image is spatially scaled 
(through an interpolation procedure) to the planetary 
image according to the wavelength ratio between the 
wavelengths of planetary and reference image. 

- The scaled reference image is subtracted from the 
planetary one. 

- Double differences: 

- Two reference images are assigned to each plane- 
tary image, with the wavelength respectively shorter 
and longer than the planetary one. The two refer- 
ence images are chosen in such a way that their 
wavelength separations from the image containing 
the planet signal are the same. 

- For each group of three images, the reference im- 
ages are spatially scaled to the planetary image ac- 
cording to the wavelength ratio between the wave- 
lengths of planetary and reference images. 

- The three images are combined acc ording to the 
double-difference formula defined by iMarois et al. I 
(120001) . 

3. The procedure at step 2 should eliminate most of the speckle 
pattern. If the pairs are selected so that the planet image is 
only present in one of the two images, the planet will not be 
canceled out. 

4. A weighted average of the cleaned differential images pro- 
vides the best final result suitable for the planet search. We 
adopted a weight for each single differential image, which is 
the reciprocal of the wavelength difference between the two 
(or three) images subtracted to obtain the considered one. In 
this way we give a greater weight to differences between im- 
ages with a smaller wavelength separation, where the speckle 
pattern has a stronger correlation. Since the planetary images 
are not scaled, the planet position will not shift with wave- 
length. 

There are three critical issues in this procedure: 

1 . To properly work we have to make an assumption about the 
spectra of the companion objects that we are looking for. 
Moreover, this works much better for spectra with large ab- 
sorption bands such as for methane-dominated planets. 

2. Each interpolation introduces noise. In our approach, the 
number of interpolations is effectively reduced to only one 
per pair (two for each group of three images when using the 
double differential imaging). 



3. Pairing of monochromatic images, and optimal weighting 
should be given according to the main noise source: 

- If errors are dominated by photon noise, the best proce- 
dure is to assign the same weights to all pairs. In this 
case, pairs should be selected to have similar (or even 
constant) wavelength separations. 

- If errors are dominated by calibration errors (speckle 
residuals), the best procedure in single differential imag- 
ing is to create pairs having the smallest possible wave- 
length separation, compatible with the gradients present 
in the planetary spectra. In this case weights should be 
assigned according to the inverse of the square of wave- 
length separation. 

- For what concerns double differences, this last approach 
is limited by the intrinsic width of the emission peaks 
in the planetary spectrum. Practically, we expect a very 
small advantage by creating groups of three images with 
the smallest possible wavelength differences. It should 
then be more advantageous to have various groups of 
three images with the same wavelength difference and 
give the same weight to all of them. 

4.2. Spectral deconvolution 

This m ethod was proposed for the first time in [Sparks & Ford I 
(l2002h and further developed in lThatte et aTl (120071) . It exploits 
that speckles are expected to change regularly with wavelength. 
Outside a given separation, defined as the bifurcation radius 
(BR), the speckle spatial excursion over the spectral range is 
larger than the planet size so that the speckle pattern associated 
to the star can be reconstructed and eliminated using regions un- 
affected by the planet image. Differently from the MDI described 
in the previous section, no assumption about the spectra of the 
companion objects is needed. 

Spectral deconvolution should offer some advantage over the dif- 
ferential imaging approach, at least outside the BR, because it 
uses the companion spectrum as a whole. The value of the the 
SPHERE IFS BR is around 0.20 arcsec for the Y-J-mode and 
about 0. 12 arcsec for the Y-H-mode. The procedure we followed 
is composed of four steps: 

- We scaled single images provided by the CSP data extrac- 
tion algorithm to a reference wavelength (in this case we 
chose the central wavelength among those of the 33 or 
38 monochromatic images). Because of this rescaling, the 
planet will be in different positions in every image. 

- We plotted the spectrum for every spaxel of the rescaled 
datacube (see Figure [3} and calculated a polynomial fitting 
function using 1/A as independent variable. The polynomial 
degree depends on the distance from the center of the im- 
age, in units of the BR. The value of this fitting function is 
then subtracted from every spectrum. The fit allows the mod- 
ulation of a given stellar halo speckle brightness with wave- 
length to be taken into account but its degree is small enough 
not to fit a potential planetary signal. This should eliminate, 
or at least reduce, the speckles or diffraction residuals. 

- The subtracted images are then rescaled back to the original 
scale according to their wavelength, in order to mantain fixed 
the planet position in all of them. 

- To search for planet signal, the three-dimensional datacube 
is collapsed to a bi-dimensional image given by the cross- 
correlation of the spectra in each spaxel with a template 
planet spectrum. This procedure enhances the signal-to noise 
of the final image. In general, in our simulations we use 
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Fig. 3. Example of the simulated spectra on a spaxel rescaled 
according to the wavelength for SPHERE IFS simulations. The 
dotted line represents the value of one single pixel on the PSFs 
used as input for the CSP code. The solid line represents the 
same thing on the images of the datacube resulting from the CSP 
simulation (passage of the light through the instrument and cor- 
rect calibration to create the datacube). The absolute values of 
these two curve have been normalized to be able to overplot one 
on the other. 



Fig. 4. Comparison between 5cr contrast obtained with single 
difference (yellow line), with multiple single differences (red 
line), and with multiple double differences (green line) from IFS 
simulations. These results were obtained for a simulation where 
no detector noise, no cross-talk, and no rotation were considered. 



which are speckles having a lifetime longer than field rotation 
but shorter than the total exposure time. 



a methane-dominated spectrum. However, as we see in 
Section [6j this procedure works well with either a flat or an 
L-type spectrum, too. 



4.3. Angular differential imaging 

In general, we assume that observations are done with the field 
fixed with respect to the IFU. In this case, the pupil rotates with 
time on an alt-az telescope, a typical value being 30° over a 
1 hour exposure time. In this frameworkQ, angular differential 
imaging (ADI) can be applied to reduce the speckle noise fur- 
ther. Vari ous codes have been written to perform ADI on real im- 
ages (see Marois et al. I l2006h . Here, we considered a variant of 
this method that we defined as azimuthal filtering ("azimuthal", 
meaning along arcs at a constant radius). This procedure is com- 
posed of the following steps: 

- For each given pixel, we searched for all spaxels that have 
similar separation (distance from the center). In our proce- 
dure, the annulus width was set at 1 pixel. 

- We plotted the value of the intensity at the selected wave- 
length for each of these spaxels against the azimuth angle. 

- We drew a fitting line through these points using a cubic 
spline curve through the averages of these points within arcs 
of length 4A/D. After various tests, we chose this value to 
avoid canceling the planet signal. 

- We subtracted the intensity value on the fitting line from the 
intensity at the selected wavelength in that spaxel. 

- The procedure was then repeated for all wavelengths. 

- The procedure was iterated over all spaxels. 

While this procedure does not completely eliminate the impact 
of static speckles, it also works well for quasi-static speckles, 

3 Angular differential imaging would actually be better applied by 
keeping the pupil static and leaving the field to rotate with respect to 
the detector. 



5. Simulations results 

In this section we review the most important results obtained 
from our simulations. As said in previous sections, we expect a 
significant improvement in the contrast using the MDI method 
compared to a simple S-SDI, when exploiting all of the many 
monochromatic images provided by an IFS. In particular we ex- 
pect that the contrast scales with the square root of the number 
of independent single differences that we can realize when us- 
ing the whole spectrum. A further improvement can be obtained 
by correctly coupling the images at different wavelengths. Since 
the contrast scales with the wavelength separation, we can pah- 
monochromatic images in ascending order of wavelength sepa- 
ration and weight them in descending order according to their 
wavelength separation. 

In Figure [4] we display these results for a simulation where no 
detector noise, no cross-talk, and no rotation were considered. 
In particular, we can see that no further gain is instead obtained 
using the multiple double-difference method. This is mainly be- 
cause realistic double differences should be made using a rather 
large wavelength separation, because of the intrinsic width of the 
methane bands. In this simulation, as in all the following ones, 
the jump in the plots around 20A/D is given by the coronagraph 
outer working angle effects. 

We can further improve the contrast obtained with our instru- 
ment by exploiting the rotation of the field with respect to the 
pupil. As seen in Figure [5] the improvement is better at large 
separations, as expected because more noise realizations can be 
sampled. If quasi-static speckles dominate, the results improve 
with the square root of the angle (and of the separation), since 
the planet images sample different noise realizations while rotat- 
ing around the stellar image. In this case the azimuthal filtering 
procedure (described above in Section PTJt can be applied. 
SD should provide better results than MDI, at least for separa- 
tions larger than the BR. This is confirmed by the plots displayed 
in Figure [6] where we show the run of the 5<x calibration limit for 
a very bright star. The case shown is for 30° field rotation with 
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Fig. 5. Run of the 5cr calibration limit with separation for three 
cases: (i) no field rotation (red line), rotation by 30° with no 
azimuthal filtering apllied (yellow line), rotation by 30° with az- 
imuthal filtering applied (green line). 
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Fig. 6. Run of the 5cr calibration limit with separation for a very 
bright star. The case shown is for 30° field rotation with az- 
imuthal filtering. Detector noises and a cross-talk with a value 
of 10~ 2 were introduced too. Red line is the result obtained with 
multiple differential imaging while the green is with the spec- 
trum deconvolution method. 



azimuthal filtering. In this case we introduced the detector noises 
using the values indicated at the end of Section [3] and a cross- 
talk total amount of 10~ 2 , too. Results obtained with the spectral 
deconvolution are slightly better than those obtained using the 
multiple differential imaging, with difference on the order of 0.2 
dex (~ 0.5 mag). As expected, better results are obtained when 
the Y-H-mode is considered. In this second case the difference 
is on the order of 0.3 dex (~ 0.7-0.8 mag), and the gain is appre- 
ciable even at a small separations (0.15 arcsec). 
From all these plots we can see that, by using SPHERE IFS and 
an appropriate combination of the analysis methods described 
above, we should be able to reach contrasts on the order of 10~ 7 
or even better at large separations from the central star. 

A synthesis of the results from our simulations is presented 
in Table[T]where we listed the contrasts obtained at different sep- 
arations from the central star using the two analysis methods for 
Y-J and Y-H-modes. 



Fig. 7. Same as Figure|6]but for the Y-H-mode. 



Table 1. 5<x calibration limit (30 degree field rotation, azimuthal 
filtering). In this table MDI stays for multiple differential imag- 
ing, while SD stays for spectral deconvolution. 









Sep. (arcsec) 




Mode 


Analysis 


0.15 


0.5 


1.0 


Y-J- 


MDI 


2.06 x 10-" 


3.09 x 10-' 


9.89 x 10- 8 


(14.21 mag) 


(16.27 mag) 


(17.51 mag) 




SD 


Inside BR 


1.64 x 10~ 7 


7.68 x 10~ 8 






(16.96 mag) 


(17.79 mag) 


Y-H- 


MDI 


1.61 x 10~ 6 


1.87 x 10~' 


1.87 x 10-' 


(14.48 mag) 


(16.82 mag) 


(16.82 mag) 




SD 


5.85 x 10- 7 


1.13 x 10~ 8 


1.55 x 10~ 7 




(15.58 mag) 


(17.37 mag) 


(17 .03 mag) 



5. 1 . Impact of Fresnel propagation 

Out-of-pupil optics could have a strong impact on the per- 
formances of any differential technique adopted in high- 
co ntrast imaging due to Fresnel propagation, as described 
by iMarois et al. I (12006). An optic that is not conjugated to a 
pupil plane will modify the light distribution in a chromatic way 
because at this location the beam intensity distribution depends 
on wavelength through diffraction effects. The closer the optic is 
to a focal plane, the larger this chromaticity. Even more severe 
is the fact that this chromaticity is no longer smooth, but cyclic 
along the spectrum, when the optic is conjugated to a height that 
is several times the Talbot length defined as 

L T = 2A 2 /A (1) 

where A is the light wavelength and A the period of a single 
sinusoidal component of the wavefront across the pupil. For an 
aberration with a given period, the pupil complex amplitude 
presenting the electromagnetic field changes from a pure wave- 
front error to a pure amplitude error over a quarter of the Tal bot 
length. Since the Talbot length is different for different periods, 
a decorrelation occurs that depends on angular separation. The 
farther an optic is from the pupil plane (in multiples of Talbot 
length), the more the decorrelation along spectral domain will 
be, and speckle correlation will be broken. 
In the case of SPHERE, the Talbot effect was expected to be 
strong for those optical components located before the lenslet 
array, such as the entran ce window, the ADC, t he derotator and 
the coronagraphic mask ( Yaitskov a et al. EoiOl) . 
To evaluate the impact of the Fresnel propagation, we cannot use 
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Table 2. Most important parameters used for the Fresnel propa- 
gation calculation. 



Opt. surface 


Con. Distance (km) 


rms WFE (nm) 


DTTS 


414 


5 


Collimator 


396 


15 


Mirror 2 


1492 


10 


Mirror 3 


4440 


10 


Field lens 


10722 


15 


0.01 


cnqu or separation 'n arcsec 
0.10 


1.00 




io"* I , , , , 

0.1 1.0 10.0 100.0 

Separation (X/D). at X= 1.31 109 

Fig. 8. Plot resulting from a simulation without the Fresnel prop- 
agation and without any rotation of the FOV (Y-H-mode). 
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Fig. 9. Same as Figure[8]but with the Fresnel propagation. 



the CAOS package that is based on the Fraunhofer pr opagation. 
We then exploited the PROPER code dKiist I l2007h to create 
new PSFs that are then used as input for the CSP code. In 
Table [2] we list the parameters used to calculate the Fresnel 
propagation in all the simulations. We report the values of the 
conjugated distance and the wavefront error (WFE) rms for all 
the considered optical surfaces. To save computing time, we 
performed all these simulations without considering the effects 
of the atmosphere (using only 1 atmospheric phase screen). Of 
course, this is not realistic, because it yields contrasts that are 
to optimistics. However, the comparison is still meaningful for 
evaluating the impact of Fresnel propagation itself. 



Table 3. Values of the contrast at different separations for sim- 
ulations with and without accounting for Fresnel propagation. 



Sep. (arcsec) 


No Fresnel 


Fresnel 


0.2 


2.40 x 10-' 


3.50 x IO"' 


0.4 


8.35 x 10~ 8 


1.31 x 10~ 7 


0.6 


1.17 x 10~ 7 


1.23 x 10~ 7 


0.8 


9.80 x 10- 8 


1.12 x 10~ 7 


1.0 


9.76 x 10~ 8 


1.19 x 10~ 7 



From comparison of the plots in Figs. |8]and|9]resulting from 
simulations that do not include and that do include the Fresnel 
propagation effects, respectively, one can see that the differences 
between the two cases are very small, as confirmed by the data 
reported in Table [3] where the values of the contrast at different 
separations for simulations performed not considering (second 
column) and considering (third column) the Fresnel propagation 
are compared. From the results of these simulations, we can then 
conclude that the effects of the Fresnel propagation are not too 
great for our instrument. This is because of the large conjugated 
distances of the optics (see Tabled and the not very large pupil 
size. Fresnel propagation is much more a concern for extremely 
large telescopes, like E-ELT (see lAntichi et al. II2010I) . 

6. Data analysis software for companion detection 
and spectral classification 

Through the simulations described in the previous sections, we 
have demonstrated the capability of the SPHERE IFS instrument 
to image extrasolar planets down to a contrast of ~ 10~ 7 at a sep- 
aration of a few tenths of arcsec. However, to fully characterize 
the newly discovered planets (temperature, chemical composi- 
tion of the atmosphere, etc.), it is important to be able to recon- 
struct their spectra at a high-fidelity level. 
To test this capability, we prepared a pipeline for a data analysis 
of the calibrated datacube resulting from our simulations. This 
procedure is composed of five different steps: 

1. Speckle noise subtraction from the original datacube using 
the spectral deconvolution algorithm. 

2. Sum of all the resulting images to create a single multi- 
wavelength image. 

3. Search for companion objects on the summed image. 

4. Extraction of a spectrum for every object found. 

5. Spectral classification of every object. 

The simulation was performed using the same PSFs as in the 
previous simulations. The FOV is rotated by 30° during the ob- 
servation and the central star is a GOV star at a distance of 10 pc 
(this corresponds to a magnitude J=3.75). Every simulation con- 
tains five planets in different positions but at the same separation 
from the central star (the planets in the same simulation are iden- 
tical). We perfomed simulations with planets at a separation of 
0.3, 0.5 and 1.0 arcsec from the central star. To avoid overlap- 
ping of the planets PSFs at the smallest separation, we decided 
to replace these single simulations with three different ones con- 
taining two different planets each (for this reason only for the 
case at 0.3 arcsec we have six planets for every single case and 
not five). We then performed different simulations with different 
luminosity contrasts between the planets and the central star, and 
adopted values of 10~ 5 , 3 x 10~ 6 , 10~ 6 , and 3 x 10~ 7 . Finally, we 
performed different simulations with different input spectra: we 
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used a late type T-dwarf spectrum (T7), an early type T-dwarf 
spectrum (T2), a late type L-dwarf spectrum (L8), and an early 
type L-dwarf spectrum (LO), taken from the spectra libraries de- 
scribed in Section loTTl To test the capability of our procedure to 
distinguish between a companion object and a background star, 
we performed simulations using a fiat spectrum (M2) at the low 
resolution of our instrument as input. Indeed, at this resolution, 
all the stellar spectra are expected to be flat. Moreover, we did 
not include any faint background galaxies because they are ex- 
pected to be spatially resolved as extended objects by our instru- 
ment. 

While all these spectra come from objects in the solar neighbor- 
hood (old objects), our results do not lose generality because, as 
shown in Section 16.2.31 the detectability of the companion ob- 
jects at a fixed effective temperature is not determined by the 
gravity effects, which in turn is the main difference between 
young and older substellar objects. 

6. 1 . Procedure description 

In this section we will describe our reduction procedure in more 
detail. The first two steps are performed using the spectral de- 
convolution method in the same way as described above in 
Section I4.2I The search for companion objects (third step) is 
composed of three different steps: 

- For each pixel of the image we compare the flux included 
in a circle centered on the analyzed pixel and the flux into 
an external annulus. The radii values of the circle and of the 
annulus can be chosen by the user, but for our analysis, we 
always adopted the values of 1.5, 2, and 4 (pixels). The user 
can choose the type of statistic to be performed on these re- 
gions: a mean or a median. From our tests we find that the 
second one was more effective in finding companion objects, 
so we always adopted it for all subsequent analysis. The pro- 
cedure finds an object if the value found for the inner circle 
is greater than for the outer annulus plus the standard devi- 
ation (on the outer annulus) multiplied by a factor that can 
be chosen by the user and that has to be considered carefully 
case by case. 

- If finding more than one object into a radius of 3 pixels, the 
procedure then retains only the most luminous. 

- Finally a two-dimensional Gaussian fit is performed on a 
small region around the newly discovered object to find its 
precise position (in 1/1000 of pixel - no evaluation of the er- 
ror on the position is done in this procedure). We try to min- 
imze the difference between the extracted PSF and the fitting 
function by performing an iterative procedure to search for 
the minimum of the difference by changing the parameter of 
the Gaussian fitting function. 

We then extracted the spectrum of the newly found object sim- 
ply summing the flux of the pixels at a distance less than 1 
pixel on every subtracted monochromatic image and subtract- 
ing from this value the median from the external annulus. We 
made the same extraction for two positions at a distance of +A/D 
from the object position along the azimuth (and then at the same 
separation of the found object) to evaluate the spectral noise. 
Subtracting the mean of these two spectra from the object spec- 
trum can then improve the final spectral classification that is the 
last step of our procedure. 

To this aim, we compared the output spectra of our simulations 
with a template spectra grid. We considered T-dwarfs from TO 
to T8 and L-dwarfs from LO to L8, with the spectral type L7 re- 
placed by L7.5, because we could not find such a spectrum in 



Table 4. Number and percentage of found objects (F.O.) and of 
spurious objects (S.O.) subdivided by the simulation input spec- 
tra for the Y-J-mode case. 



Sp. Type. 


F.O. 


% F.O. 


S.O. 


% S.O. 


T7 


55 out of 64 


85.9 


7 out of 62 


11.3 


T2 


61 out of 64 


95.3 


24 out of 85 


28.2 


L8 


61 out of 64 


95.3 


29 out of 90 


32.2 


LO 


60 out of 64 


93.8 


31 out of 91 


34.1 


M2 


59 out of 64 


92.1 


24 out of 83 


28.9 



the literature. In addition, we also considered F1V, GOV, K5V, 
M2V and M8V type star spectr a. The data fo r the T -dwarfs 
temp late spectra were taken from lLooper et al. I (120071) for TO, 
from Burgasser et al. (2004) for spectra from Tl to T5 and for 
T8 and from Burgasser at al. (200g) fo r T6 and T7. The da ta for 
the L-dwarfs spectra were taken from iTesti et al. I (l200lh . The 
stellar spectra have been taken from the IRTF online Spectral 
Library (http : //irtfweb . if a . hawaii . edu/~ spex/IRTF_ 
Spectral. Library/index. html). 

The spectral classification was obtained by a cross-correlation 
(using the IDL routine C_CORRELATE) between the output 
spectrum of each simulation and the template spectra. The spec- 
tral type with the highest cross-correlation coefficient is the one 
assigned to the simulated planet. 

6.2. Results 

6.2.1. Companion detection 

In Tables 3] and |6l we display the numbers and the percentages 
of found objects divided according to the spectral type of the 
input spectra of the simulations (second and third columns) for 
the Y-J and the Y-H modes respectively. In the fourth and in the 
fifth columns of the same tables, we instead report the numbers 
and the percentages of the spurious objects found with our pro- 
cedure. It is apparent that we are able to find most of the simu- 
lated objects both for the Y-J and the Y-H modes, but the method 
works better in the second case. Moreover, the number of spu- 
rious objects found is much lower for the Y-H-mode case than 
in the Y-J-mode case. We do stress that almost all the simulated 
objects that we are not able to find in the final image are for the 
cases at a separation of 0.3 arcsec where the background noise 
from the central star is greater. Indeed, it is 100% complete for 
companions down to a contrast of 3 x 10~ 7 and separations of 
0.5 arcsec, while it is complete at more than 90% for contrasts 
other than the worst case with a contrast 3 X 10~ 7 and separation 
0.3 arcsec. 

To confirm this we report in Tables |5] and [7] the number and 
the percentage of those found and of the spurious objects as in 
Tables [4] and [6] but, in this case divided according to the lumi- 
nosity contrast of the simulated objects. It is apparent from these 
tables that we are able to find almost all the simulated objects 
down to a contrast of 10~ 6 , while we lose more than 25% of the 
simulated objects with a 3 x 10~ 7 luminosity contrast using the 
Y-J-mode. On the other hand, we are able to find more than 90% 
of the objects with a contrast of 3 x 10~ 7 using the Y-H-mode. 

6.2.2. Spectral classification 

In Figure[l0]we show the spectral classification of all real objects 
found, and the spectral classification of the spurious objects. For 
the real objects, we can see three high peaks corresponding to the 
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Table 5. Number and percentage of found objects (F.O.) and of 
spurious objects (S.O.) subdivided by different contrasts of the 
simulated objects for the Y-J-mode case. 



Contrast 


F.O. 


% F.O. 


S.O. 


% S.O. 


10- 5 


80 out of 80 


100.0 


44 out of 124 


35.4 


3 x 1(T 6 


80 out of 80 


100.0 


26 out of 106 


24.5 


1(T 6 


78 out of 80 


97.5 


16 out of 94 


17.0 


3 x 10- 7 


58 out of 80 


72.5 


29 out of 87 


33.3 


ble 6. Same as Table [U but for the Y-H-mode. 


Sp. Type. 


F.O. 


% F.O. 


S.O. 


% S.O. 


T7 


63 out of 64 


98.4 


1 out of 64 


1.6 


T2 


64 out of 64 


100.0 


out of 64 


0.0 


L8 


61 out of 64 


95.3 


12 out of 73 


16.4 


LO 


63 out of 64 


98.4 


22 out of 85 


25.9 


M2 


61 out of 64 


95.3 


14 out of 75 


18.7 



3x1 Cr 6 
10"' 
3x1 0" 1 







F1 GO K5 M2 MS L0 L1 L2 L3 L4 L5 L6L7.5L8 TO T1 T2 T3 T4 lb T6 T7 TS 
Spectral type 

Fig. 11. Same of Figure [l0]but with the histograms divided for 
different luminosity contrast. The dashed lines represents the 
spurious objects. 



Table 7. Same as Table[3] but for the Y-H-mode case. 



Contrast 


F.O. 


% F.O. 


S.O. 


% S.O. 


io- 5 


80 out of 80 


100.0 


25 out of 105 


23.8 


3 x 10- 6 


79 out of 80 


98.7 


14 out of 93 


15.1 


io- 6 


78 out of 80 


97.5 


6 out of 84 


7.1 


3 x 10- 7 


75 out of 80 


93.7 


4 out of 79 


5.1 



F1 GO K5 M2 M8 L0 L1 L2 L3 L4 L5 L6L7.5LS TO T1 T2 T3 T4 T5 T6 T7 T8 
Spectral type 

Fig. 10. Histogram with the number of objects (red) and of spuri- 
ous objects (blue) found for every spectral type in the Y-J-mode 
case. 



M8, Tl, and T8 spectral types. The M8 peak is given by the con- 
tribution of objects with both an M2 and an L0 input spectrum. 
The Tl peak is given by the objects with L8 and T2 input spec- 
tra. In this case, however, the peak is quite low and the objects 
classification is more dispersed. Finally, the T8 peak is given by 
T7 input spectra objects. In general, apart from the case of the 
L0 spectral type, it seems that our procedure tends to classify the 
objects with later spectral types rather than the actual ones. We 
do not have any particular peak in the final distribution for the 
spurious objects. 

A similar histogram, but for the Y-H-mode, is displayed in 
Fig. [12] Even in this case, we have three peaks in the distribu- 
tion of the real objects. The first one again corresponds to the 



Fl GO K5 M2 M8 L0 L1 L2 L3 L4 L5 L6L7.5L8 TO T1 T2 T3 T4 T5 T6 T7 T8 
Spectral type 

Fig. 12. Same as Figure[10l but for the Y-H-mode case. 



3x1 0" s 

10"' 

3x10"' 



F1 GO K5 M2 US L0 L1 L2 L3 L4 L5 L6L7.5LB TO T1 T2 T3 T4 lb 16 17 TS 
Spectral type 

Fig. 13. Same as FigureQT] but for the Y-H-mode case. 



M8 spectral type and comes from the contribution of the simu- 
lations with M2 and L0 input spectra objects. This means that, 
as for the Y-J-mode case, these two spectral types seem to origi- 
nate in a degeneracy. The second peak is around the T4 spectral 
type and is mainly given by the T2 input spectra simulations, 
but from the L8 simulations too. The L8 simulations do not give 
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Table 8. Cross-correlation coefficients considering the effects of 
the gravity with Y-J-mode. 







log(g) = 4.0 


log(g) = 5.5 


log(g) 


= 4.0 


0.88 


0.40 


log(g) 


= 5.5 


0.80 


0.51 



Table 9. Cross-correlation coefficients considering the effects of 
the gravity with Y-H-mode. 







log(g) = 4.0 


log(g) = 5.5 


log(g) 


= 4.0 


0.89 


0.77 


log(g) 


= 5.5 


0.72 


0.90 



in general a correct identification. Indeed, these objects are rec- 
ognized alternatively as L2 type or early T type. The last peak 
is at the T7 spectral type and it is given exclusively by the T7 
simulations objects (the T8 detections are given by the simula- 
tions with a separation of 0.3 arcsec). In Figs.QT|and[T3]we dis- 
play the same histograms, divided in these cases according to the 
different luminosity contrasts of the simulated objects indicated 
with different colors in the figures. From these figures we can 
see that the overall distribution of the spectral classification is 
very similar to the global one displayed in Figs. [T0land[T2ldown 
to a contrast of 10~ 6 , while they are much more dispersed for the 
3 x 1CT 7 contrast where the spectral classification becomes much 
less effective. 

6.2.3. Effects of the gravity 

To further test the capability of our procedure to distinguish dif- 
ferent objects, we performed different simulations using as in- 
put the synthetic spectrum of one object with T e ff = 
and log(g) = 4.0 and of another one with the same tempera- 
ture and log(g) = 5.5. All the simulations were performed for 
five different objects (with the same characteristics) at a separa- 
tion from the central star of 0.5 arcsec and a contrast of 3 x 10~ 6 . 
Furthermore, we performed simulations both for the Y-J and the 
Y-H-modes. 

For the simulations with the Y-J-mode, all the objects with 
log(g) = 4.0 were recognized as T8 spectral type (with values 
of the cross-correlation coefficients around 0.75) while the ob- 
jects with log(g) = 5.5 were recognized as T7 (4 cases) or T8 
(1 case). In this second case, the values of the cross-correlation 
coefficients are on the order of 0.77. 

On the other hand, for the simulations with the Y-H-mode all 
the objects with log(g) = 4.0 were recognized as T8 spectral 
type but with higher values of the cross-correlation coefficients 
(more than 0.93), while all the objects with log(g) = 5.5 were 
recognized as T6 spectral type (cross-correlation coefficients on 
the order of 0.92). 

In Tables [8] and [9] we report the values of the mean coefficients 
from the cross-correlation between the output and the input spec- 
tra for the Y-J-mode and for the Y-H-mode, respectively. From 
these results it is apparent that, in the case of the Y-H-mode we 
are able to correctly classify the objects for the gravity effects, 
while for the Y-J-mode all the simulated objects are classified as 
logfe) = 4.0, . 

In conclusion, from our analysis it seems that the Y-H-mode 
is the best solution for correctly distinguishing between objects 
with different gravities. 



7. Conclusions 

We performed detailed simulations of the performances of the 
SPHERE IFS instrument and considered different data analysis 
methods that can be exploited to reduce data coming from the 
instrument. In particular, we exploited the multiple spectral dif- 
ferential imaging (MDI), the spectral deconvolution (SD), and 
the angular differential imaging (ADI). This latter seems to be 
especially useful associated with one of the other two methods. 
It turned out that SD is slightly more effective in reducing the 
speckle noise than the MDI, and it is less sensitive to the charac- 
teristics of the planetary spectrum. From our analysis, however, 
we can now conclude that, in the best cases, the IFS channel 
of SPHERE should be able to image companion objects around 
nearby stars down to a contrast of almost 10~ 7 at a few tenths of 
an arcsec. 

We then performed detailed simulations to test the possible im- 
pact of Fresnel propagation on the final performances of the in- 
strument. This issue created some concerns especially the pres- 
ence of optics before the lenslet array but, from our simulations 
made under the same IFS optical setup, a negligible difference 
results between the achievable constrasts with or without con- 
sidering the effects of Fresnel propagation. 
Because the SD method, as said above, allows better results, we 
used it to perform a new analysis fn the capability of the instru- 
ment to find and to characterize companion objects of the central 
star. 

We then prepared a pipeline with the aim of reducing the dat- 
acube resulting from our simulations. To test the effectiveness 
of this procedure in finding and characterizing planets, we per- 
formed a series of simulations with different companion objects' 
input spectra, different separations, and different contrasts be- 
tween the simulated planets and the central star. From these sim- 
ulations and exploiting the spectral deconvolution method com- 
bined with some ADI, we were able to image extrasolar planets 
down to a luminosity contrast with respect to the central star of 
3 X 10~ 7 . In this way we confirmed the results obtained with the 
previous run of simulations. 

We have generally been able to find almost all the simulated 
objects at the larger separation considered (0.5 and 1.0 arcsec), 
while the method is less effective at a separation of 0.3 arcsec. 
However, even in this case, we were able to find more than the 
90% of the simulated objects using the Y-J-mode and more than 
the 95% of the objects using the Y-H-mode. 
For the spectral reproducibility of our procedure, we can adopt 
the following conclusions: 

- The greater the separation from the central star, the larger the 
possibility to reconstruct the planets' spectra with precision 
(considering planets with the same luminosity contrast). 

- Planets with greater luminosity contrast more easily have a 
precise spectrum reconstruction. 

- This method allows us to reconstruct and to classify the T 
type spectra very well while spectral reconstruction and clas- 
sification seem to be less precise for earlier spectral types. 
However, even in these cases, the spectral classification gen- 
erally has a precision of a few spectral types (4 or 5 in the 
worst cases). 

- Stellar spectra (M and earlier spectral types) are clearly 
distinguished from T-type spectra while some ambiguity is 
present for L-type companions. This implies that, in most 
cases, the characterization of the nature of the detected ob- 
jects (companion vs. field star) can be obtained from discov- 
ery data alone without waiting for common proper motion 
confirmation. Ambigous cases of L-type companions can be 
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disentangled in several cases from the properties of the ob- 
ject and the parent star (e.g. L-type companion are expected 
only above a given contrast threshold). 

- The Y-H-mode allows a better spectral classification than for 
the Y-J-mode. 

- For what concerns the effects of the gravity, they are better 
disentangled using the Y-H-mode than using the Y-J-mode. 
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